;+
; NAME:
;       GENERATE_LPSTR
;
; PURPOSE:
;       Genertate the lpstr structure (LOOP STRUCTURE)
;
; CATEGORY:
;       IRDC distance modeling subroutine
;
; CALLING SEQUENCE:
;       GENERATE_LPSTR, ROTC=rotc, SUFF=suff
;
; INPUTS:
;       NONE
;
; OPTIONAL INPUTS:
;       ROTC -- 0 = Reid, 1 = IAU, 2 = CLEM
;
; KEYWORD PARAMETERS:
;       NONE
;
; OUTPUTS:
;       NONE
;
; OPTIONAL OUTPUTS:
;       SUFF -- Filename suffix to identify rotation curve used,
;               needed by some routines
;       RTIT -- Plot title, needed by some routines
;
; COMMON BLOCKS:
;       FFORE_BLOCK
;
; MODIFICATION HISTORY:
;
;       Created:  05/02/12, TPEB -- Initial version, finally got sick
;                                   of pasting code block into various
;                                   routines.
;
;-

PRO GENERATE_LPSTR, ROTC=rotc, SUFF=suff, RTIT=rtit
  
  IF n_elements(rotc) EQ 0 THEN rotc=0
  
  ;; Define the COMMON BLOCK
  COMMON FFORE_BLOCK,n,rb3,R0,d,R,Z,tau,f_data,sig_f,corr,farlist,do_Tfit,$
     lpstr,rho_hi,rho_h2,rho_star
  FORWARD_FUNCTION BGPS_FFORE_MODEL, BGPS_FFORE_LOOPIE, STAR_DENSITY
  
  
  RESOLVE_ROUTINE, 'BGPS_FFORE', /COMPILE_FULL_FILE, /EITHER, /NO_RECOMPILE
  
  n = n_elements(rb3)
  
  ;; Get galactic params
  defsysv, '!MW', exists = exists
  IF NOT exists THEN galactic_params 
  ;; Using !MW parameters -- Goes into COMMON block!
  R0 = !MW.R0
  d = dindgen(!MW.NBINS*5.)*!MW.BINSIZE + !MW.BINSTART
  
  message,'Setting the stage (re-used calculations)...',/inf
  ;; Do the calculations needed by every loop and place them into a
  ;; structure to save on computation time
  lpstr = REPLICATE( CREATE_STRUCT('kdn',0.,'kdf',0.,$
                                   'r',dblarr(!MW.NBINS*5.),$
                                   'z',dblarr(!MW.NBINS*5.),$
                                   'rho_hi',dblarr(!MW.NBINS*5.),$
                                   'rho_h2',dblarr(!MW.NBINS*5.),$
                                   'rho_star',dblarr(!MW.NBINS*5.)), n)
  
  ;; Do distances for the particular curve
  CASE rotc OF
     0: BEGIN                             ;; Reid et al. (2009)
        FOR i=0L, n-1 DO BEGIN
           
           lpstr[i].kdn = KDIST(rb3[i].l, rb3[i].b, rb3[i].vlsr, /NEAR)
           lpstr[i].kdf = KDIST(rb3[i].l, rb3[i].b, rb3[i].vlsr, /FAR) 
           
           omni_lbd2rz, rb3[i].l, rb3[i].b, d, r, z
           lpstr[i].R = r
           lpstr[i].Z = z
        ENDFOR
        suff = ''
        rtit = 'Reid et al. (2009)'
     ENDCASE
     1: BEGIN                             ;; IAU standard
        FOR i=0L, n-1 DO BEGIN
           
           lpstr[i].kdn = KDIST_IAU(rb3[i].l, rb3[i].b, rb3[i].vlsr, /NEAR)
           lpstr[i].kdf = KDIST_IAU(rb3[i].l, rb3[i].b, rb3[i].vlsr, /FAR) 
           
           R0 = 8500.d
           omni_lbd2rz, rb3[i].l, rb3[i].b, d, r, z, R0=R0
           lpstr[i].R = r
           lpstr[i].Z = z
        ENDFOR
        suff = '_iau'
        rtit = ''
     ENDCASE
     2: BEGIN                             ;; Clemens (1985)
        FOR i=0L, n-1 DO BEGIN
           
           lpstr[i].kdn = KDIST_CLEM(rb3[i].l, rb3[i].b, $
                                     rb3[i].vlsr, /NEAR)
           lpstr[i].kdf = KDIST_CLEM(rb3[i].l, rb3[i].b, $
                                     rb3[i].vlsr, /FAR) 
           
           R0 = 8500.d
           omni_lbd2rz, rb3[i].l, rb3[i].b, d, r, z, R0=R0
           lpstr[i].R = r
           lpstr[i].Z = z
        ENDFOR
        suff = '_clem'
        rtit = 'Clemens (1985)'
     ENDCASE
  ENDCASE
  
  ;; Calculate the HI & H2 densities along the LOS for each source
  ;;   here, and place into the lpstr for later retrieval (speed-up).
  FOR i=0L, n-1 DO BEGIN
     lpstr[i].rho_hi = hi_density(lpstr[i].R,lpstr[i].Z)
     lpstr[i].rho_h2 = h2_density(lpstr[i].R,lpstr[i].Z)
     ;; The 2400 = 2 * (H1 + H2) is the normalization
     lpstr[i].rho_star = star_density(lpstr[i].R,lpstr[i].Z)
     
  ENDFOR
  
END
